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Abstract 

The variability of AVHRR calibration coefficients 
in time was analysed using algorithms of linear and 
non-linear time series analysis. Specifically we have 
used the spline trend modeling, autoregressive pro- 
cess analysis, incremental neural network learning 
algorithm and redundancy functional testing. The 
analysis performed on available AVHRR data sets 
revealed that (1) the calibration data have non- 
linear dependencies, (2) the calibration data depend 
strongly on the target temperature, (3) both cal- 
ibration coefficients and the temperature time se- 
ries can be modeled, in the first approximation, as 
autonomous dynamical systems, (4) the high fre- 
quency residuals of the analysed data sets can be 
best modeled as an autoregressive process of the 
10th degree. 

We have dealt with a non-linear identification prob- 


lem and the problem of noise filtering (data smooth- 
ing). The system identification and filtering are sig- 
nificant problems for AVHRR data sets. The al- 
gorithms oulined in this study can be used for the 
future EOS missions. Prediction and smoothing al- 
gorithms for time series jf calibration data provide 
a functional characteiization of the data. Those al- 
gorithms can be particularly useful when calibration 
data are incomplete or sparse. 

1. Introduction 

EOS-Earth Observing System 
The suite of instruments to be flown on the Earth 
Observing System is intended to provide a compre- 
hensive data set of global observations of the Earth 
in a broad range of sensor wavelenghts. Global cov- 
erage sensors will be the primary tools for data col- 
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lecting. Those instruments will play a major role 
in the global studies of our environment under the 
EOS project. To provide a quality data set for sci- 
entists, the raw data sets need to be calibrated. Ob- 
viously this effort requires a good calibration algo- 
rithm. We will focus in this paper on analysis of 
time series of calibration coefficients. This analysis 
can help in defining an optimum value of calibra- 
tion coefficient based on the previous and current 
values determined from on-board calibration mea- 
surements. Also, the description of calibration data 
sets in terms of e.g. nonlinear dynamical systems 
and autoregressive processes will help in improving 
both the short and longterm estimation accuracy of 
calibration coefficients. 

AVHRR Calibration Coefficients Data 

Data sets from Advanced Very High Resolution 
Radiometer (AVHRR) sets were analysed. The 
AVHRR system was flown on the National Oceano- 
graphic and Atmospheric Administration (NOAA) 
operational meteorological satellites (NOAA-9/11). 
At a given time, the calibration data from one chan- 
nel consists of a pair of numbers, the slope and inter- 
cept, that represent a straight line, the calibration 
curve. For one channel we obtain two time series, 
the slope as a function of time, and the intercept as 
a function of time. We have analysed the available 
data sets; year 1988, day 25; and year 1990, day 
183. 

The Calibration Time Series 
The calibration data sets form time series. A time 
series {y(<)} can be thought of as a realization of a 
stochastic process. A stochastic process can be de- 
scribed as a sequence of random variables. The best 
studied are linear stochastic processes. We used 
standard modeling techniques [Box, Jenkins 1976], 
based on autoregressive (AR) models to analyse the 
high frequency components of the data. Among our 
long term goals is to establish a short term, medium 
term and long term model for the estimate of cal- 
ibration coefficients for AVHRR data. We would 
like to approximate calibration time series by linear 


models whenever we can avoid using computation- 
ally costly non-linear models. Testing for nonlinear- 
ity is not a trivial task. Recently tests for nonlinear- 
ity based on the redundancy (entropy) functionals 
have been introduced [ Palus, 1993 ] as a tool for de- 
tecting chaotic dynamical systems. We used redun- 
dancy functionals to detect nonlinear dependencies 
of the time series in our data sets. Knowing that 
we are dealing with non-linear models we would like 
to build a nonlinear prediction model for calibration 
coefficients time series as a multivariable function of 
external variables (the principal one is undoubtedly 
the target temperature). At this phase of the anal- 
ysis we have used an incremental (recursive) neural 
network architecture to simulate the dynamics of 
the behavior of calibration coefficients and the tar- 
get temperature for one orbit worth of data. 

We have observed four main signal components in 
the data sets. The first component, a slow trend, 
corresponds to the aging of the sensors. The second 
trend, a nonlinear ( pseudo-periodic ) dependence, 
corresponds mainly to the effect of the day or night 
part of the spacecraft orbit on the sensors ( the sen- 
sors are housed in the spacecraft which is affected 
by variable external pseudo-periodic conditions ) 

The third component consists of the medium fre- 
quency signal (a few minute period). 

The calibration coefficients are contaminated by 
noise. 

2. The Time Series Modeling 
Goals 

Data Averaging 

The first goal of our analysis of the AVHRR data 
sets was to define the current value of a calibration 
coefficient, i.e. the value at a given time. Since 
a data set of calibration coefficients is a time se- 
ries contaminated with noise, this task consists of a 
trend estimation, by a trend modeling and a data 
smoothing procedure. The results of this process 
are time series of calibration coefficients with noise 


filtered out .This process includes the estimates of 
confidence intervals. 

Analysis of Residuals 

The second goal, and a part of the data averaging 
procedure, is analysis of the residuals. As the series 
can have autocorrelated terms, it needs to be mod- 
eled by an autoregressive process. 

Calibration Data Dynamics 

The third goal was to investigate the modeling and 
predictibility of the calibration coefficients in time 
in terms of a nonlinear dynamical system 

dx/ dt — f(x(t),t, u(<)) (1) 

where x(t) represents a time series, /() is an 
unnknown non-linear function which we want to 
model and approximate, u(t) represent external pa- 
rameters (e.g. the target temperature). 

Dynamics Diagnostic 

A number of tests for non-linearity have been pro- 
posed [Tong 1993]. We have used a redundancy- 
based test for non-linearity [Palus 1993]. 

AVHRR Applications 

Applications of modeling techniques for the AVHRR 
data sets are far reaching. Firstly, they will allow 
us, on a rigorous basis, to define the current cali- 
bration values of sensors more accurately and sta- 
tistically fully qualified in terms of confidence in- 
tervals. Secondly, they will improve the long term 
and short term estimates of calibration coefficients 
through capturing the system dynamics. 

3. Calibration Value Estima- 
tion Procedure 

For AVHRR data sets, we made the assumption that 
the observed series can be described by the model 


x(t) = q(t) + v(t) + e(t) (2) 

where t is the discrete reference time, q(t) repre- 
sents a long-term trend (pseudo-oscilations with pe- 
riod equal approximately to 102 minutes, or 1 orbit), 
v(tf) is a component describing the short-term peri- 
odicities. Finally, e(t) is white noise. The first term 
q(t) is driven by external phenomena (e.g. tempera- 
ture, light intensity ), and the dynamics of the sen- 
sors. 

A number of methods are available for signal compo- 
nent modeling. We focused in the present study on 
the method of splines, autoregressive modeling, re- 
dundancy functionals and adaptive neural network 
methods. 

Trend Modeling in the Time Domain 
(Splines) 

The first task, before we can analyse noise, consists 
of estimating and removing the component q(t) from 
the equation (2). A splines smoothing is the stan- 
dard modeling choice [Wegman 1983]. As the char- 
acter of the component q(t) changes several times 
during one period, an algorithm is needed to cap- 
ture this change in the trend. The number of knots 
can be controlled by the Akaike’s AIC criterion [Eu- 
bank, Speckman 1990]. 

Dynamics Modeling in the Time Domain 
(Dynamical Systems) 

The physical mechanisms governing the evolution 
of the system of calibration coefficients are not fully 
understood but we assumed that the data is pro- 
duced by an underlying generator which is a low 
dimensional dynamical system. 

The first component q(t) is driven by external phe- 
nomena (e.g. temperature, the light intensity). This 
relation is complex and non-linear. In the first ap- 
proximation we assumed that the 2nd order au- 
tonomous difference equations 

dxjdt = /(*(<)) (3) 

can simulate a dynamical system that governs 
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the calibration coefficients. In other words we as- 
sumed that the calibration coefficients are fully de- 
scribed by the dynamics of the temperature evolu- 
tion, which is , for simplicity, described by a 2nd 
order autonomous syatem. The 2nd order calibra- 
tion coefficient system generates a nonlinear map 

g : S C R? — ► R 
x(i+ 1) = <7(x(z),x(i- 1)) 

where the set S contains the training exemplars. 
This map can be captured by a neural network al- 
gorithm. 

Dynamics Modeling in the Time Domain 
(Neural Networks) 

We have simulated one orbit worth of the 1988 data 
set by the Cascor learning algorithm. Because we 
assumed the 2nd order dynamics (see the previous 
section) we set the number of inputs for our neu- 
ral net to two. Tha trained network creates a delay 
map characterizing the underlying dynamical sys- 
tem. Many different architectures of feedforward 
neural networks have been used for signal model- 
ing and prediction [Chen, Billings, 1992]. Nonlin- 
ear system identification and modeling using neural 
networks has become very popular tool for mod- 
eling and identification of non-linear autonomous 
systems [Vemuri, 1994]. The neural network algo- 
rithms offer the flexibility of infinitely many parame- 
ters. The number of parameters is usually indirectly 
controlled through cross-validation. The method of 
cross-validation consists, in practical terms, of two 
phases. First a net is trained on a subset of data, 
secondly the trained net is tested on another data 
set. This process is repeated until the residual er- 
ror on the test data hits its minimum. In our nu- 
merical simulation the cross validation test has been 
implemented using the Cascor code [Fahlman,1992]. 
The follow-on tests will demonstrate the numerical 
agreement of the prediction and the real data. The 
real data for additinal orbits are being prepared for 
neural network prediction tests. 


Series Prediction by Bootstrapping 

The trained net can be used for predicting the future 
values of the time series. We have used trained nets 
to predict various dynamic systems. Many authors 
have used different neural network architectures to 
predict chaotic systems. The standard test cases are 
logistic equation, Mackey-Glass nonlinear differen- 
tial delayed equation and van der Pol equation. 

The methodology of predicting, known as boot- 
strapping, works as follows: For an n input, m out- 
put network, a training exemplar is formed by tak- 
ing n+m consecutive values from the data series to 
be extrapolated (predicted). Starting at an arbi- 
trary value x(£), the first n values ( x(i ), x(n-\- i — 
1)) are presented to the network inputs. The target 
values are the next m values, thus a general training 
exemplar for an n input, m output network can be 
represented as (z(i), x(n + i— 1)..., x{n -f m + i- 
1)). Each successive exemplar is formed by start- 
ing one value beyond the previous starting value. 
The number of exemplars sufficient to provide the 
desired accuracy is to be determined by numerical 
experiments. After the iterative learning procedure 
has converged we have a map (e.g. for one output, 
m=l) 

x(i + 1) = g(x(i ), .... x(i - n + 1)) 

By bootstrapping the net into the future , this 
map can be iterated to give 

x(i+2) = j 7 (flf((xi), ...,x(i-n+l)),x(i), ...,x(i-n+2)) 

The trained network can be analysed and differ- 
ence equations describing the underlying dynamics 
can be recovered [ Lowe, Webb 1994 ]. 

Noise Analysis in the Time Domain (AR- 
Autoregression) 

The AR is the standard method for modeling lin- 
ear dynamical systems. In the next step we have 
to decide about the character of the second compo- 
nent v(t). The standard methodology recommends 
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analysing spectrum and periodogram, in order to 
decide between the trigonometric regression or au- 
toregression. The FFT transfomation, Fig. 3, did 
not show well isolated high-frequency components 
of our signal. Therefore the standard methods us- 
ing low-pass filters for noise removal are not directly 
applicable for our data sets. The periodogram re- 
vealed several isolated peaks in some regions. This 
feature indicated that trigonometric functions might 
be a good candidate for modeling of our data set. 
However, we have not observed this behavior of pe- 
riodograms through the entire time domain. On the 
other hand, the partial autocorrelations showed sig- 
nificant dependencies of the lagged vector compo- 
nents for lags up to 10, less significant dependencies 
for the higher lags. That is why the AR approach 
has been preferred for statistical analysis of the se- 
ries. In other words, if we denote the residuals of 
the first step of analysis r(t) = x(f) — q(t ), we have 
the model 

jfc 

r(t) = v(t ) + e(t) = ^2 b(i)r(t - i) + e(t), 

i = l 

where b(i) are the parameters and k is the order 
of the autoregressive process AR(k). 

Nonlinearity Test 

The real-world data show usually some degree of 
nonlinearity. Our test was based on information- 
theoretic (redundancy) functionals. The redun- 
dancy test is based on the fact that noise, linear and 
nonlinear structures of a time series are represented 
by qualitatively different redundancy functionals. 

Testing for nonlinearity of the calibration coeffi- 
cients time series was performed using redundancies. 
Those methods were recently proposed for testing 
of non-linearity of dynamical systems [Palus, 1992] 
and are based on the general concepts of informa- 
tion theory. 

The linear redundancy L(Xi y ... t X n ) of an n- 
dimensional random variable with zero mean and 


covariance matrix C is defined as 

n n 

L(Xi, ...,X n ) = l/2£>g( Cii ) - l/2^1og( ( r I ) 

i= 1 i=l 

where c,, are the diagonal elements and cr t are the 
eigenvalues of the n x n covariance matrix C . 

The n-dimensional (non-linear) redundancy is de- 
fined as 

R{X u ... y X n ) = H(Xi) + ....H(X n ) + H(X Xl „.X n ) 

These two functionals, as the functions of the time 
lag r, provide measures that differentiate linear and 
non-linear structures presented in the lagged ver- 
sions of the component x(i). The lagged version 

(z(f), x(t + r), ...x(t + (n - l)r)) 

of x(t) is a realization of the random variable 
(X\, ". } X n ) where n is the embbeding dimension. 

4. Data Characterization 

Visual inspection of the data from different time 
periods (1988,1990) and for different platforms re- 
vealed common features. We have searched the 
available data sets for features that could explain 
the variability of the calibration coefficients. Two 
features were outstanding. The long-term compo- 
nent u(/) behaves, in the first approximation, as 
an autonomous dynamical system within one orbit, 
Fig. 6. 

Also PRT counts (Platinum Resistance Thermome- 
ter), located in the word 20 of the 103- word of 
HRPT minor frame [Kidwell. 1991] behave as an au- 
tonomous dynamical system. The target tempera- 
ture is temperature of the internal target. This tem- 
perature can be calculated from the output of four 
PRT counts located in words 18, 19, 20 of HRPT 
minor frame. The conversion of PRT counts c to 
absolute temperature is accomplished by 

T(K) = J2 a jc i 

j= o 
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. The patch temperature ( word 21) has been con- 
stant for the inspected data sets. 

5. Summary of Numerical re- 
sults 

Our numerical results and observations are summa- 
rized in figures 1 through 17. In the following the 
first two digits indicate the year, and the next three 
digits Julian day. The slope coefficients vary signif- 
icantly over one orbit. Fig.l shows the slope coef- 
ficient series for approximately one period (88025). 
Fig. 2 shows the PRT counts series for 88025. The 
FFT for the slope data (88025) did not show any sig- 
nificant isolated frequencies, in Fig. 3. Fig. 4 shows 
the typical relation between a slope series and PRT 
counts for the same time period (88025). The PRT 
dynamical system delay diagram (88025) shown in 
Fig. 5 is a typical example of a map which can be 
learned by a neural network. The slope dynamical 
system delay diagram (88025) is another example of 
the dynamics learnable by a net (Fig. 6). The PRT 
series (90183) shows similar features to the series 
(88025) even when the time difference between these 
series is several years (Fig. 7). The same can be said 
about the slope values (90183), in Fig. 8. Fig. 9 is 
a slope vs. PRT diagram. A data smoothing proce- 
dure is illustrated in the following three figures. Fig. 
10 shows a series of intercept coefficients (88025). 
Fig. 11 shows the smoothed intercept (88025) data. 
The intercept residuals (88025) showed in Fig. 12 
were modeled by an autoregressive process. Fig. 13 
shows the best fit of the 10th order AR(10) (88025). 
Tests for non-linearity are captured in the remain- 
ing four figures: Fig. 14 shows slope linear redun- 
dancy statistics, dim=2 (90183). Fig. 15 shows 
slope linear redundancy statistics, dim=3 (90183), 
Fig. 16 shows slope non-linear redundancy statis- 
tics, dim=2 (9018315), Fig. 17 shows slope non- 
linear redundancy statistics, dim=3 (9018315). The 
original slope data sets were first differentiated. The 
results for different lags are very similar. In both 


cases the difference between the linear and nonliner 
measures is quite obvious (a factor of 10). That 
clearly indicates a need for a nonlinear description 
of the AVHRR data sets. 

6. Conclusions 

We have theorized, based on our numerical simula- 
tions, that the dynamics of the coefficients (within 
a time frame of days) is almost a periodic process 
with random input variables, due to such random 
effects as variable cloud cover. 

A possible future application target of the pro- 
posed algorithms is short term, medium term and 
long term data prediction for MODIS. The MODIS 
instrument is scheduled to be launched in 1998. One 
of the important tasks in the processing of MODIS 
data will be to determine the most accurate value of 
the calibration coefficients and their corresponding 
uncertainties, in other words to model and identify 
time series generated by MODIS sensors in different 
wavebands. This can be accomplished by combining 
the predited values furnished by a nonlinear (neural 
net, nonlinear AR) or linear model. A smoothing 
algorithm provides calibration values without noise. 
A prediction model will be useful especially when 
we dealing with the sparse or missing data. An- 
other potential application is the detection of sud- 
den degradation , as opposed to gradual aging, of 
a sensor. That will be characterized by an abrupt 
change in the dynamics of the sensor. To detect 
this change we may compare predicted values of the 
time series, generated by the sensor, with the new 
observed values. A discrepancy between the pre- 
dicted and observed values flags a faulty sensor. 


7. Future Work 

We are in the process of using the proposed method- 
ology for more extensive data sets. The more ex- 
tensive data sets will allow us to demonstrate the 



predictive power of the algorithms and to establish 
their error bounds. Special attention will be given 
to the integration algorithms for data from differ- 
ent calibrators ( Solar Diffuser, Spectral Radiomet- 
ric Calibration Assembly, Black Body, Space View), 
see chapter 5 in [ Guenther, 1994 ]. Theoretical re- 
search will cover experimenting with different neu- 
ral network architectures (cascade, RBF, incremen- 
tal architectures). A special effort will be devoted 
to recovering non-linear difference equations from 
the trained neural network and to establishing error 
bounds of nonlinear predictors. 
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